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Abstract. We calculate the ionisation fraction in protostellar disk models using a number of different chemical 
reaction networks, including gas-phase and gas-grain reaction schemes. The disk models we consider are conven- 
tional a-disks, which include viscous heating and radiative cooling. The primary source of ionisation is assumed 
to be X-ray irradiation from the central star. For most calculations we adopt a specific disk model (with accretion 
rate M — 10 -7 Moyr -1 and a = 10 -2 ), and examine the predictions made by the chemical networks concerning 
the ionisation fraction, magnetic Reynolds number, and spatial extent of magnetically active regions. This is to 
aid comparison between the different chemical models. 

We consider a number of gas-phase chemical networks. The simplest is the five species model proposed by 
Oppenheimer & Dalgarno (1974). We construct more complex models by extracting species and reactions from 
the UMIST data base. In general we find that the simple models predict higher fractional ionisation levels and 
more extensive active zones than the more complex models. When heavy metal atoms are included the simple 
models predict that the disk is magnetically active throughout. The complex models predict that extensive re- 
gions of the disk remain magnetically uncoupled ("dead") even when the fractional abundance of magnesium 
iMg = 10~ 8 . This is because of the large number of molecular ions that are formed, which continue to dominate 
the recombination with free electrons in the presence of magnesium. 

The addition of submicron sized grains with a concentration of x gT = 10 -12 causes the size of the "dead zone" 
to increase dramatically for all kinetic models considered, as the grains are highly efficient at sweeping up the 
free electrons. We find that the simple and complex gas-grain reaction schemes agree on the size and structure of 
the resulting "dead zone", as the grains play a dominant role in determining the ionisation fraction. We examine 
the effects of depleting the concentration of small grains as a crude means of modeling the growth of grains 
during planet formation. We find that a depletion factor of 10 -4 causes the gas-grain chemistry to converge to 
the gas-phase chemistry when heavy metals are absent. When magnesium is included a depletion factor of 10 -8 
is required to reproduce the gas-phase ionisation fraction. This suggests that efficient grain growth and settling 
will be required in protoplanetary disks, before a substantial fraction of the disk mass in the planet forming 
zone between 1-10 AU becomes magnetically active and turbulent. Only after this has occurred can gas-phase 
chemical models be used to predict reliably the ionisation degree in protoplanetary disks. 

Key words, accretion, accretion disks - MHD - planetary systems: protoplanetary disks - stars: pre-main sequence 



1. Introduction 

Observations of young stars in star forming regions point 
to the ubiquity of protostellar disks. These disks show 
evidence for active accretion, with the canonical mass 
flow rate being 10 -8 M Q yr _1 (e.g. Sicilia-Aguilar et al. 
2004) , requiring a means of transporting angular momen- 
tum in the disks. At the present time only one mechanism 
has been shown to work: MHD turbulence generated by 
the magnetorotational instability (Balbus & Hawley 1991; 
Hawley & Balbus 1991; Hawley, Gammie & Balbus 1996). 

There are questions, however, about the applicability 
of the MRI to cool and dense protostellar disks, as the ion- 
isation fraction is expected to be low (e.g. Blaes & Balbus 



1994; Gammie 1996). Magnetohydrodynamic simulations 
of disks including ohmic resistivity (Fleming, Stone & 
Hawley 2000) show that for magnetic Reynolds numbers 
Re m smaller than a critical value -Re™' , turbulence cannot 
be sustained and the disks return to a near-laminar state. 
The value of i?e mt depends on the magnetic field configu- 
ration. For net-flux vertical fields turbulence is sustained 
when Re m > 100. For zero-net flux fields, turbulence dies 
when Re m < 10 4 . Although there remains some uncer- 
tainty about the precise value of -Re™' due to the possible 
role of additional non-ideal MHD effects, recent simula- 
tions by Sano & Stone (2002) indicate that the inclusion 
of Hall E.M.F.s has little impact on the value of -Re™ 1 *. A 
value of -Re™ 1 * = 100 typically corresponds to a gas-phase 
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electron fraction of x[e~] ~ 10 -12 for most disk models. 

The question of whether protostellar disks can sustain 
MHD turbulence has important implications beyond the 
issue of how mass accretes onto the central protostar. In 
particular turbulence may play an important role in planet 
formation. The growth of planetesimals is controlled by 
the velocity dispersion (e.g. Wetherill & Stewart 1993), as 
is the growth of protoplanetary cores (Kokubo & Ida 2000; 
Thommes et al. 2003). Simulations of planets in disks in- 
dicate that turbulence may have a significant influence 
on planet migration and gap formation, as well as the 
velocity dispersion of planetesimals (Nelson & Papaloizou 
2003, 2004; Winters, Balbus & Hawley 2003; Nelson 2005). 
Obtaining detailed knowledge of the structure of turbulent 
protoplanetary disks is a crucial part of developing models 
of planet formation. 

A number of studies of the ionisation fraction in pro- 
tostellar disks have appeared in the literature. The early 
work by Gammie (1996) put forward the idea that disks 
may have magnetically "active zones" sustained by ther- 
mal or cosmic ray ionisation, adjoining regions that were 
magnetically inactive - "dead zones" . More recent stud- 
ies have examined this issue, but in each case different 
chemical reaction networks and/or disk models have been 
employed, making comparison between them difficult. For 
example, Gammie (1996) used a layered a-disk model, 
and considered cosmic-rays as the primary ionising source. 
Sano et al. (2000) used the Hayashi (1981) minimum mass 
solar nebula model in conjunction with a more complex 
chemical model that included dust grains. Glassgold et 
al. (1997) and Igea et al. (1999) used the minimum mass 
model and considered X-rays as the most important ion- 
isation source. Fromang et al. (2002) considered the ion- 
isation fraction in conventional a-disk models, and fo- 
cussed on the role of heavy metals in a simple reaction 
network. Matsumura & Pudritz (2003) adopted the exter- 
nally heated, passive disk model proposed by Chiang & 
Goldreich (1997) in conjunction with the Sano et al. (2000) 
chemical reaction network. Recent work by Semenov et al. 
(2004) studied chemistry in the D'Alessio et al. (1999) disk 
model, and employed a complex gas-grain chemical net- 
work including most reaction from the UMIST data base. 

In this work we study the chemistry and ionisation 
fraction in conventional a-disk models, similar to the ap- 
proach taken by Fromang et al. (2002). We include X-ray 
ionisation from the central star as the primary (external) 
ionisation source, and neglect cosmic-ray ionisation. One 
of our primary goals is to compare the ionisation fraction 
predicted by the different reaction networks, within the 
framework of the same underlying disk model. 

We study both gas-phase chemical models, and gas- 
grain models. We begin by examining the ionisation frac- 
tion generated by the simple gas-phase, five species model 
proposed by Oppenheimer & Dalgarno (1974) for study- 
ing the electron fraction in dark clouds. We study also 
a number of more complex gas-phase chemical networks. 
These include the Sano et al. (2000) kinetic model (ap- 
plied to the gas phase), and networks constructed from 



various species sets that are drawn from the UMIST data 
base (Le Teuff et al. 2000). In general we find disagreement 
between these models, with the more complex networks in 
particular producing lower fractional ionisation than the 
simpler schemes. This is because a large number of molec- 
ular ions form that are able to recombine quickly with free 
electrons. Even a relatively large abundance of heavy met- 
als (magnesium) is unable to render the disks fully active 
in these cases. 

We construct a series of gas-grain chemical models 
that are generalisations of the gas-phase networks. In line 
with expectations, we find that the addition of submicron 
sized (r gr — 10~ 5 cm) grains causes the magnetically ac- 
tive zone to shrink dramatically, with each of the models 
being in basic agreement about its size and structure. As a 
simple approach to modelling the effects of growth during 
planet formation, we examine the effects of depleting small 
grains on the chemistry. We find that depletion factors of 
10 -4 are required to reproduce the gas-phase chemistry 
for models in which heavy metals are absent. Depletion 
factors of ~ 10~ 8 , however, are required to reproduce 
the ionisation fraction obtained when metals are included. 
The effectiveness of very small numbers of grains in modi- 
fying the equilibrium ionisation fraction depends crucially 
on there being a constant supply of neutral grains, these 
being generated through recombination reactions between 
negatively charged grains and positive ions. These models 
suggest that efficient growth of small grains will be re- 
quired before they can safely be neglected from chemical 
networks that predict the ionisation fraction in protoplan- 
etary disks. 

This paper is organised as follows. The a-disk mod- 
els we construct are described in section The chemical 
models are described in section [31 and our scheme for cal- 
culating the X-ray ionisation rate is described in sectionQJ 
The results of the calculations are presented in sectional 
and extensions to the basic models are presented in sec- 
tion Finally we summarise the paper in sectional 

2. Disk models 

Numerous studies of the ionisation structure in protoplan- 
etary disks have appeared in the literature, with a vari- 
ety of disk models being used. The resulting ionisation 
structure differs from model to model, as they differ in 
their density and temperature profiles. The early study 
by Gammie (1996) used a layered a-disk model. Sano et 
al. (2000) used the Hayashi (1981) minimum mass solar 
nebula model, as did Glassgold et al. (1997) and Igea et 
al. (1999). Matsumura & Pudritz (2003) adopted the ex- 
ternally heated, passive disk model proposed by Chiang 
& Goldreich (1997), and recent work by Semenov et al. 
(2004) studied chemistry in the D'Alessio et al. (1999) 
disk model which includes U.V. heating from the central 
star and internal viscous heating. In this work we study 
the chemistry and ionisation fraction in conventional a— 
disk models, similar to the approach taken by Fromang et 
al. (2002). 
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Fig. 1. This contour plot shows the density distribution 
for our a = 1(T 2 , M = 1CT 7 Moyr" 1 disk model. The 
contour lines refer to values of 10 -11 and 10 -10 gem -3 , 
respectively. 

A detailed description of the disk model and mod- 
elling procedure is given in Ilgner et al. (2004). Here we 
just sketch out the important details. We compute a disk 
model with inner and outer radii [0.1, 10] AU orbiting 
with Keplerian velocity about a central star with mass 
1 M Q . The vertical structure is obtained at a given ra- 
dius by solving for hydrostatic and thermal equilibrium. 
Radiative and convective transport is included with opac- 
ities taken from Semenov et al. (2003), and disk heating 
occurs through viscous dissipation alone. With these in- 
gredients and appropriate boundary conditions the disk 
model is uniquely specified for given values of a and the 
mass accretion rate M. Most of our models are computed 
with a — 10~ 2 and M — 10~ 7 M yr~ 1 , giving a mass 
of 0.0087 M Q within the computational domain. Contour 
plots showing the density and temperature distributions 
are shown in figures ^ and [21 We consider a model with 
a = 5 x 10~ 3 and M = 10~ 8 M yr _1 toward the end of 
this paper. Similar models and procedures for calculating 
a-models were used by Papaloizou & Terquem (1999) and 
Fromang et al. (2002). 

3. Chemical models 

We use standard techniques to solve the kinetic equations 
that describe the chemical models we consider in this pa- 
per. We aim to provide sufficient detail to enable an inter- 
ested reader to reproduce our results, so full descriptions 
of our procedures are presented in appendices. In this sec- 
tion we begin by describing important aspects of the ki- 
netic models we have implemented, prior to discussing rate 
coefficients and initial abundances. 

3.1. Kinetic models 

A kinetic model is defined by the components (or species) , 
the chemical reactions that cause time dependent changes 
in the abundances of species, and by the underlying ki- 
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Fig. 2. This contour plot shows the temperature distribu- 
tion for our a = 10~ 2 , M = 10~ 7 Moyr" 1 disk model. 
The contour lines refer to values of 100, 500, and 1000 K, 
respectively. 



netic equation for each component. The kinetic equations 
are given by a set of stiff ordinary differential equations 
which are solved numerically. We have confirmed both al- 
gebraically and numerically that conservation of elements 
and charge is satisfied by our kinetic models. We have im- 
plemented the following reaction schemes, and examined 
the resulting ionisation fractions in a-disk models. 

3.1.1. Oppenheimer & Dalgarno model 

Oppenheimer & Dalgarno (1974) introduced a simplified 
reaction scheme to approximate the equilibrium electron 
abundance in dense interstellar clouds. Using an analytic 
approximation, this scheme has been used to calculate the 
electron fraction in protoplanetary disks by a number of 
authors (Gammie 1996; Glassgold et al. 1997; Fromang et 
al. 2002). 

For a given total particle density n, Oppenheimer & 
Dalgarno introduced a simple kinetic model involving two 
elements ("m" and "M"), five components ("m", "m + ", 
"M" , "M + " , and "e~" ), and four chemical reactions. Here, 
m and m + represents a molecule and its ionized counter- 
part, while M and M + refer to a heavy neutral metal atom 
and its ionized counterpart. e _ denotes the free electron 
in the gas phase. 

In this model, ions are produced by ionisation and a 
charge-transfer reaction, while dissociative recombination 
and radiative recombination are assumed to destroy ions. 
The reaction set is represented schematically in table ^ 
where £, a, /3, and 7 are the ionizing flux, the dissociative 
recombination rate coefficient, the charge-transfer rate co- 
efficient, and the radiative recombination rate coefficient. 

One can derive the equilibrium electron abundance an- 
alytically using a simple approximation (which is n ~ 
n [m] where n[m] denotes the number density of the neu- 
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Table 1. List of reactions (and the corresponding rate co- 
efficients) considered by Oppenheimcr & Dalgarno (1974). 
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tral molecule m). The equilibrium electron concentration 
£oo[e - ] is given by the roots of the algebraic equation 1 
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where x[M] is the equilibrium particle concentration of 
metals. Neglecting the contributions due to the heavy 
metal atom equation (JTJ can be reduced to its quadratic 
form whose solution is 
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In the opposite limit, i.e. assuming a:oo[M] 3> Xoo[c _ ], 
similar expressions can be derived depending on further 
assumption being made. See for example equation (27) 
in Oppenheimer & Dalgarno (1974) and equation (14) in 
Fromang et al. (2002). 

3.1.2. Umebayashi & Nakano model 

Sano et al. (2000) calculated the ionisation fraction in 
a minimum mass solar nebula disk model adopting a 
reaction scheme which was investigated originally by 
Umebayashi & Nakano (1990) and Nishi et al. (1991) un- 
der dense interstellar cloud conditions. 

This reaction scheme extends the more schematic gas- 
phase scheme of Oppenheimer & Dalgarno (1974) by intro- 
ducing additional neutral and ionized gas-phase species, 
and gas-grain interactions. Sano et al. (2000) considered 7 
ions (H + , ELj", , He + , C + , m + , M + ) and grain particles 
with seven different grain charges (gr, gi^, gr 2± , gr 3± ). 
The ions m + and M + represent a molecular and metal 
ion, respectively. 

The reactions considered by Sano et al. (2000) are 
listed in table|2which is taken from their paper. Note that 
the ions 0+, 0+, OH+, 2 H+, CO+, CH+, and HCO+ are 
represented by the molecular ion m + (Sano, private com- 
munication). Our only omission was triple charged grains, 
which Sano et al. (2000) included. 

The results obtained by Sano et al. (2000) can be 
used for testing our implementation of this kinetic model. 
Applying the same physical conditions (which includes a 
surface density distribution S(i?) based on the minimum 
mass solar nebula model), the same ionisation flux, and 
rate coefficients (including the value of S c - — 0.6 for all 



For simplicity we dropped the symbol oo from x[e ] and 
x[M]. 



Table 2. List of reactions considered by Sano et al. (2000) 
excluding triple charged grains. Notation: X + denotes the 
ions of the kinetic model: [H + , H^~, H^", He + , C + , m + , 
M + ] while Y, Y a< j s refer to neutral gas-phase species and 
adsorbed species onto grain particles, respectively. 
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temperatures and grain charges, where 5 C - denotes the 
sticking probability of electrons onto neutral or charged 
grains), and the same equilibrium values of the abun- 
dances of neutral components 2 as Sano et al. (2000) used, 
we are able to reproduce their results. The equilibrium dis- 
tributions of the most important components are shown in 
figure[3]and they correspond extremely well to the fiducial 
model presented in figure (1) of Sano et al. (2000). This 
indicates that our implementation of the Umebayashi & 
Nakano reaction scheme is correct. 
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Fig. 3. Equilibrium abundances (per hydrogen nucleus) of 
some representative particles at the disk midplane, which 
should be compared to the fiducial model of Sano et al. 
(2000). The distribution of e" and of HCO+ is given by 
the solid and the dashed line, respectively. 

3.1.3. Extended kinetic models 

The kinetic models introduced by Oppenheimer & 
Dalgarno and by Umebayashi & Nakano may serve as pro- 
totypes for gas-phase chemistry and gas-grain chemistry, 
respectively. Based on these prototypes we have developed 
more complex kinetic models by extending the number of 
elements, species, and reactions. 

We have formalised our approach to developing ex- 
tended kinetic models as follows. We start with all gas- 
phase species listed in table IA.l1 This species set contains 
the elements H, He, C, O, N, S, Si, Mg and Fe. We ex- 
tract reactions involving these species from the UMIST 
database (Le Teuff et al. 2000). The same species set was 
used by Willacy et al. (1998), Willacy & Langer (2000), 
and Markwick et al. (2002) in their studies of disk chem- 
istry. This forms the basis of the gas-phase chemical net- 
work that we implement. 

We also consider chemical evolution due to gas-grain 
collisions and grain-grain collisions. As yet there is no 
UMIST-like database that describes these types of reac- 
tions, so we have developed our own model. We use a 
nomenclature that differentiates between reactions that 
involve mantle species and those that do not. We refer to 
the former as "mantle chemistry" and the latter as "grain 
chemistry". We now describe the mantle chemistry, fol- 
lowed by the grain chemistry. 

Consider first the neutral components of the previously 
described gas-phase chemistry. We define a new class of 
species: mantle species - which are the adsorbed counter- 
parts of the gas-phase neutral species, with the exception 
of Helium which is assumed to remain in the gas phase. 
These reactions are represented schematically in table [3J 
where "X" is a neutral species and the "gr" are grain par- 
ticles whose charges are denoted by their indices. A colli- 
sion between a neutral species X and a grain with charge q 
simply results in the formation of a mantle species X(gr 9 ) 



Table 3. List of generalised reactions determining the 
mantle chemistry. Neutral gaseous species are represented 
by the symbol X, its ionised counterpart by X + . The nota- 
tion X(gr) and X(gr~) refer to neutral gas-phase species 
adsorbed onto grain particles with no excess charge and 
with single negative excess charge, respectively. 
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which physically represents atoms/molecules of species X 
sitting on grains with charge q. 

We now consider the "grain chemistry" . The reactions 
are represented schematically in table^J Reactions 1 and 2 
show that a positive ion colliding with a negative grain re- 
sults in the ion being neutralised, but remaining in the gas 
phase, and the grain being left with one less charge - pro- 
vided a neutral counterpart X of the ion X + exists in the 
UMIST database. Collisions between electrons and grains 
neutralise positive grains, and negatively charge neutral or 
singly-negative charged grains. The model incorporates 
single and double charging of grains. Grain-grain colli- 
sions result in conservation of grain number but involve 
charge exchange. There are some species in the set de- 
scribed in tablc lA.il whose collisions with grains cannot be 
described by these reactions. In particular these are colli- 
sions between grains and those ions X + (such as Hjj" and 
NH^) with no neutral counterpart. The products of these 
collisions are unknown. By neglecting these particular col- 
lisions between ions and grains one may overestimate the 
ionisation fraction. To remedy this we have assumed that 
collisions between ions X + and negatively charged grains 
form the same products as the dissociative reactions in 
the gas phase: 

X+ + gr" — > Y + Z + gr 
X+ + e- — > Y + Z 

Quite a few of the dissociative recombination reactions in 
the UMIST database are branch reactions. We did not 
use the same branching ratio for the corresponding gas- 
grain collision. Instead, we assumed that each branch is 
equally weighted and that the weights are normalised, i.e. 

We note that gas-grain models have been considered 
by Willacy et al. (1998), Willacy & Langer (2000) and 
Markwick et al. (2002). These authors, however, did not 
include charged grains in their models. The models com- 
puted by Sano et al. (2002) did include charged grains. 

We have computed a number of models using gas- 
phase chemistry only, and gas-phase plus mantle and 
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Table 4. List of generalised reactions determining the 
grain chemistry. X + denotes those components of the 
species set which have a neutral counterpart X. 
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Table 5. List of kinetic models. Models with indices 1 and 
5 refer to both prototypes discussed in the text. The re- 
maining models refer to the extended models. While mod- 
els with indices lower than 4 refer to pure gas-phase chem- 
istry only, models with indices 4, 6, and 7 include grain 
and mantle chemistry. 
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grain chemistry. We now describe the implementation of 
these different models. 

3.1.4. Model setup 

The kinetic models we have constructed and applied are 
listed in table |3J Note that for the fractional abundance 
of grains x gY — > 0, model4 and model7 reduce to modell 
and model3, respectively. 

We considered three different kinetic models of the 
gas-phase chemistry: 

1) modell : The first and the simplest one is the reaction 
scheme introduced by Oppenheimer & Dalgarno (1974). 
Species M can be represented by any neutral metal atom, 
and species m represents molecular hydrogen. The model 
is described in detail in subsection 3.1.1. 

2) mode!2 : This kinetic model used two elements, simi- 
lar to the Oppenheimer & Dalgarno model. All reactions 
involving hydrogen and magnesium only were extracted 
from the UMIST database. Note that the reactions in- 
volved here differ from the Oppenheimer & Dalgarno 
model. 

3) mode!3 : All species listed in table IA.1I were included, 
giving the most complex kinetic model we have studied 
so far. The model includes the elements: H, He, C, O, 
N, S, Si, Mg, and Fe. The same set of gaseous species was 
adopted by Willacy et al. (1998), Willacy & Langer (2000) 
and Markwick et al. (2002) in their studies of disk chem- 
istry. 

In addition, we considered four different kinetic models 
which combine gas-phase chemistry with gas-grain pro- 
cesses: 

1) model4 : The simplest one is a modification of 
Oppenheimer & Dalgarno's model (modell) obtained by 
adding the grain and mantle chemistry. The elements are 
(molecular) hydrogen and magnesium, which need to be 
specified to define the evaporation temperature from grain 
surfaces. 

2) mode!5 : This refers to the kinetic model by Umebayashi 
& Nakano (1990) as applied by Sano et al. (2000). It is de- 
scribed in subsection 3.1.2 . 



3) mode!6 : This model was designed to improve the ki- 
netic model of Sano et al. (2000) by including all species 
listed in table lA. il eontaining the 5 elements (H, He, C, O, 
and Mg). We extracted 79 gaseous species from table IA.l1 
which include free electrons. 31 of the 78 gaseous species 
have an adsorbed neutral counterpart on the grains, so 
making up the mantle species. Evolving the 79 gaseous 
species requires the extraction of 799 gas-phase reactions 
from the UMIST database. Applying our schemes for grain 
and mantle chemistry (see tables |3] and 0J) 291 additional 
reactions were introduced. 

4) mode!7 : Including all the species listed in table IA.1I 
we construct the most complex kinetic model of gas- 
phase, grain, and mantle chemistry we have studied so far. 
Following the same procedure which was used to set up 
model6, we have evolved 174 gaseous species (including 
free electrons) and 69 mantle species. The same species 
set (gaseous + mantle species) was used by Willacy et 
al. (1998), Willacy & Langer (2000) and Markwick et al. 
(2002). 1965 reactions were extracted from the UMIST 
database accompanied by 621 additional reactions involv- 
ing grains. 

3.2. Rate coefficients 

Expressions for the rate coefficients can be obtained from 
macroscopic considerations. The mean collision rate coef- 
ficient is given by < er c v > where cr c and v denote the 
cross section and the velocity of the molecules. While the 
rate coefficients of the gas-phase reactions can be calcu- 
lated easily using information in the UMIST database (Le 
Teuff et al. 2000), computation of the rate coefficients for 
gas-grain processes requires more detailed information. 
We discuss this further in the following sections, but note 
here that we assume a grain radius of r gr = 10~ 5 cm in 
all calculation described in this paper. 
Compared with previous versions of the UMIST database, 
its latest version (Le Teuff et al. 2000) accounts for the 
temperature dependence and specifies the valid temper- 
ature range for all reactions. Where explicit information 
about the temperature range is not available, the rate co- 
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Fig. 4. Sticking probability of an electron onto neutral 
grains and charged grains as a function of the gas- 
tcmpcrature T. Curves are labeled by values of v — Ze/q. 
D and r denotes the polarization potential depth D of the 
grain and the grain radius, respectively, while Ze and q 
refer to the grain charge and the charge of the gas-phase 
particle. 

efficients have been attributed a temperature range of 10 
- 300 K. Since the parameterisation relies on a mathe- 
matical fit to complex data over a restricted temperature 
range, the validity of the rate coefficients must be checked 
when used outside their temperature range. We used con- 
servative estimates when calculating the (temperature- 
dependent) rate coefficients: if the rate coefficient de- 
creases with temperature T, we always used the low tem- 
perature limit T\ as a cut-off value; for lower tempera- 
tures the rate coefficient was approximated by T = T\. 
If the rate coefficient increases with temperature, we al- 
ways used the high temperature limit T2 as a cut-off value; 
for higher temperatures the rate coeffficient was approxi- 
mated by T = T2. 

3.2.1. Collisions between ions/electrons and grains 

The rate coefficients were estimated following the paper 
by Draine & Sutin (1987). Compared to the approxima- 
tion of Umebayashi & Nakano (1980) which is used widely 
in astrochemistry, Draine & Sutin included contributions 
from the Coulomb potential and the polarization of the 
grain particles induced by the Coulomb field. 

After some algebra a dimensionless reduced rate coef- 
ficient J can be extracted from the general expression for 
the rate coefficient. Depending on the ratio v of the charge 
Ze of the grain particle to the charge q of the gas-phase 
particle, which relates to the strength of attractive and 
repulsive Coulomb interactions, an approximate formula 
for the reduced rate coefficient can be derived. We applied 
the equations given by Draine & Sutin (1987), especially 
equation (3.3) for v — 0, equation (3.4) for v < 0, and 
equation (3.5) for v > 0. 

When calculating the rate of gas-grain particle inter- 
action, we took account of the fact that not every colli- 



Fig. 5. Sticking probability of neutral species onto grains 
as a function of the gas-temperature T. AE S denotes the 
amount of energy transferred to the grain particle while 
D is the dissociation energy. The other values of D/k-Q are 
1000, 1500, and 5000 K. 

sion leads to sticking to the grain surface by introducing 
an average probability S - the "sticking coefficient" . We 
calculated the sticking coefficients - which depend on the 
energy, the mass, the charge of the incident gas-phase par- 
ticle, the grain charge and the radius of the grain particle - 
by applying equations (B5), (B6), (B13) which are given in 
Nishi et al. (1991). While the electron sticking coefficient 
S c - depends strongly on the temperature and strength of 
the Coulomb interaction (see figure BJ, the sticking coef- 
ficients S*x+ for ions depend only weakly on the tempera- 
ture. For simplicity, we assumed Sx+ = 1 which is a good 
approximation. 

3.2.2. Collisions between neutral gas-phase particles 
and grain particles 

We have assumed (see reactions listed in table |3J) that 
collisions between neutral gas-phase particles and grain 
particles form mantle species. In order to specify the rate 
coefficients for these collisions, neutral gas-phase particles 
are considered to stick onto grain surfaces due to thermal 
adsorption only. The rate coefficient can be obtained by 
averaging the kinetic energy Ei of the incident gas-phase 
particles over the thermal energy distribution, assuming 
that the cross section a c is independent of Ei and is given 
by the geometrical cross section only. The rate coefficient 
due to thermal adsorption is then k — a c < v > where 
< v > denotes the mean thermal velocity. 

When calculating the sticking rate we accounted for 
the fact that more than one collision is needed to trap the 
gas-phase particle onto the grain surface. The procedure 
for calculating the sticking probability <Sx of neutral gas- 
phase species is similar to that for calculating S c -. The 
probability P e that a particle with kinetic energy Ei will 
be adsorbed is approximated by P e — exp{— e 2 /2}, where 
e is given by e = Ei / ^JDAE S (see Hollenbach & Salpeter 
1970). D and AE S denote the dissociation energy and the 
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amount of energy transferred to the grain particle due to 
lattice vibration, respectively. For each neutral gas-phase 
component we approximated D by its binding energy for 
physical adsorption Ed (see table [A~2"|) : the value of A_E S 
was approximated by 2.0 x 10 -3 eV. 
Figure |S] serves as an illustration for the temperature de- 
pendency of the sticking probability Sx for different values 
of dissociation energy. 

3.2.3. Collisions between grain particles 

The rate coefficients for the four different types of grain- 
grain collisions (see table 0} were calculated by applying 
equation (3) of Umebayashi & Nakano (1990). Here the 
value of the mass density of the grain material p gI is taken 
to be pg r = 3 g cm~ 3 . 

3.2.4. Desorption processes 

Desorption processes cause the ejection of mantle species 
from the grain particles. Here, we considered the thermal 
desorption of mantle species only. The average time for a 
mantle species to overcome the surface binding is given 
by equation (2) in Hasegawa et al. (1992), replacing the 
potential energy barrier between adjacent surface poten- 
tial energy wells, Eh, by the binding energy for physical 
adsorption Ed- The binding energy for physical adsorp- 
tion Ed for each considered mantle species is listed in 
table IA.21 The standard value of 1.5 x 10 15 cm 2 for the 
surface density of sites is assumed. 

3.3. Initial abundances 

While the composition of the gas may vary depending on 
the specific kinetic model, the total particle density n of 
the gas phase can be considered as a conserved quantity. 
The variation in the mean molecular weight p of the gas 
due to chemical reactions is negligible as the abundances 
of molecular hydrogen and helium are approximately con- 
stant. Hence the value p — 2.33 is kept constant (assuming 
a solar composition of the most abundant elements H and 
He), and the total particle density of the gas phase can 
be approximated by n — p/(pmn). Here p and toh denote 
the mass density of the gas and the mass of the hydrogen 
nuclei, respectively. 

We have adopted the same elemental abundances used 
by previous authors in their studies of the ionisation frac- 
tion (Sano et al. 2000, Fromang et al. 2002, Semenov et 
al. 2004). The abundances of those elements which are not 
covered by the model of Sano et al. (2000) are taken from 
Semenov et al. (2004). Hydrogen and helium excepted, 
all elements are assumed to be depleted in the gas phase 
compared with solar abundances (shown in table [Bj. The 
fraction 8 of the elements in the gas phase is: 1 for H and 
He, 0.2 for C and O, 2.2 x 10" 1 for N, 2.7 x 10~ 4 for Si, 
and 4.9 x 10~ 3 for S. While the reservoir of all the other 
elements is fixed by these values, the abundances of met- 



Table 6. Elemental abundances per hydrogen nucleus. 



Element 


Abundances 


H 


1.00 


He 


9.75 x 10~ 2 


c 


3.26 x 10~ 4 


N 


1.12 x 10~ 4 


O 


8.53 x 10~ 4 


Mg 


variable 


Si 


3.58 x 10~ 5 


S 


1.85 x 10~ 5 


Fe 


variable 


gr 


variable 



als (Mg and Fe) are taken as free parameters and will be 
specified elsewhere. 

The correct choice of initial abundances of chemical 
species is uncertain as we have no detailed knowledge of 
how material in molecular clouds evolves as it collapses to 
form a protostellar disk, or how it evolves during the early 
stages of disk formation and evolution. We used the sim- 
plest approach for tackling this issue. With exception of 
hydrogen (which is completely locked in molecular form), 
the elements are taken to be atomic. The number of hydro- 
gen nuclei nn (which is needed to get the species' particle 
concentrations) is approximated by nu = ^/2~3i c i> where 
Ci refers to the elemental abundances (per hydrogen nu- 
cleus) of element i in the gas phase. Note the difference 
in the notation: xpQ] denotes the particle concentration 
of species normalised to the total particle density n. 
Xi denotes the elemental particle concentration of the el- 
ement i. 

Since we consider grain particles with different electric 
charge excess, only neutral grains were assumed to exist 
initially, i.e. the total number density n gr of grain particles 
is given by n gr = n[gv\ with nfgr^ = n[gr 2± ] = 0. 

4. Sources of free electrons 

In this work we treat the free electrons e _ as an inde- 
pendent component of the reaction network and solve a 
kinetic equation for their abundance a:[e _ ]. 

Many elementary ionisation processes are considered 
in the UMIST database. Reactions which may cause the 
formation of free electrons are: i) photoreactions ii) reac- 
tions due to cosmic-ray particles which includes photore- 
actions induced by cosmic-ray particles and iii) chemiion- 
isation reactions. 

We have neglected the effects of the stellar and inter- 
stellar ultraviolet radiation field, which drive the photore- 
actions. There is great uncertainty about the amount of 
unattenuated UV flux that reaches the disk surface. In ad- 
dition the exponential decay of this field with increasing 
column density means that it can only act as source of 
free electrons in a thin skin at the disk surface. As we are 
primarily interested in calculating the depth of the mag- 
netically active zones in the disk models, this thin surface 
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Table 7. 



1. 


H 2 


— H+ 


+ 


e 




0.97 


Ceff 


2. 


H 2 


— » H+ 


+ 


H 


+ e" 


0.03 


Ccff 


3. 


He 


— » He+ 


+ 


e~ 




0.84 


Ccff 



layer is unimportant. 

Galactic cosmic-ray particles are unlikely to reach the 
surface of the inner regions of a protostellar disk due 
to their exclusion by T Tauri winds. Indeed the solar 
wind blocks cosmic-ray particles with energies less than 
100 MeV from entering the inner solar system. Preferring 
a conservative estimate of the electron abundance, we did 
not consider any of the 11 cosmic-ray particle reactions 
listed in the UMIST database. For the same reason, we 
did not consider any of the cosmic-ray-induced photore- 
actions. 

Only one chcmiionisation reaction, i.e. 

CH + O — ► HCO+ + e~ 

is tabulated in the UMIST database; the importance of 
this particular type of (associative) ionisation for the 
formation of carbon monoxide in interstellar clouds was 
pointed out by Dalgarno et al. (1973). 

Other additional sources of free electrons such as the con- 
tributions due X-ray photons, thermal ionisation, and the 
decay of radioactive elements are analysed in the subse- 
quent sections 4.2. 5.1, and 6.0. 

4.1. X-rays 

We have included the effects of X-ray radiation on 
the chemistry. Since the first X-ray observatories were 
launched, the X-ray properties of YSOs have been well 
documented (see the review by Glassgold et al. 2000). 

We adopted values Lx = 10 30 erg s _1 and fceTx = 
3 keV for the total X-ray luminosity Lx and the plasma 
temperature Tx, respectively, corresponding to a young 
star of approximately solar mass. Adopting the scheme 
described by Fromang et al. (2002), we modelled the X- 
ray source as a coronal ring of radius 10 R© centered on 
the star. 

4.2. X-ray ionisation rate 

We assumed that photoelectric absorption of X-ray pho- 
tons dominates while neglecting Compton scattering. For 
the X-ray photon energies considered, this assumption 
simplifies the description of the interaction of X-ray pho- 
tons with matter in a way that depends on the attenuat- 
ing column densities. In addition, the results presented by 
Igea & Glassgold (1999) indicate that the difference in the 
total X-ray ionisation rates between the two cases "are 
about the same at small and large values of the vertical 
column density and differ by about one order of magni- 
tude at the scattering shoulder" . 




0.1 1 5 10 

r [AU] 



Fig. 6. The effective X-ray ionisation rate Ceff per hydro- 
gen nucleus (including contributions due to thermal ioni- 
sation of potassium). The disk parameters are a = 10 -2 
and M = 10~ 7 M yr _1 . The contour lines refer to values 
of Coff: 1CT 19 , 10~ 21 , and 10~ 23 s~ x . 

Following Krolik & Kallman (1983), we modelled the 
X-ray flux of the optically thin bremsstrahlung spectrum 
by the exponential law F e = Fo exp{— e/fceTx}- For larger 
X-ray optical depths 3 tx, the local X-ray flux F per 
unit energy is attenuated by F = F e exp{— Tx}- We fur- 
ther assumed that the X-ray ionisation only depends on 
the elemental composition of the gas and is independent 
of the species composition. Hence, the absorption of X- 
ray photons can be simplified by introducing an effec- 
tive photoionisation cross section (per hydrogen nucleus) 
o'eff = Tlii^iCTi. Here Oi refers to the partial photoionisa- 
tion cross section of each atomic species and Xi denotes 
the particle concentration of the element i. We used the 
expression for the effective (total) photoionisation cross 
section given by Igea & Glassgold (1999) for the gas-phase 
depletion of the heavy elements compared to solar abun- 
dances. 

We calculated the effective X-ray ionisation rate Ceff 
for hydrogen by assuming that Ceff is dominated by sec- 
ondary ionisation. The number of secondary ionisations 
N sec per unit of primary photoelectron energy is consid- 
ered to be independent of energy. Finally, the effective 
X-ray ionisation rate (per hydrogen nucleus) Ceff is ap- 
proximated by 

/>oc 

Ccff « N scc / F(E)a cS dE . (3) 

The value of N sec is about 26 for hydrogen if the energy E 
is measured in units of keV. £7 m i n denotes the minimum 
X-ray photon energy. Hence, for values of photon energies 
less then E m [ n the photon will be absorbed before the disk 

3 The X-ray optical depth is given by rx = a(E)N. a is the 
total gas X-ray phot oabsorpt ion cross section per hydrogen 
nucleus while N denotes the (particle) column density along 
the line of sight toward the X-ray source. 
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surface is reached. Our choice of E nlin = 0.1 keV - which 
is above the required energy for ionisation of H2 (15.4 eV) 
- accounts for absorption in the source (see Shang et al. 
2002). 

The integral can be transformed into a dimensionless 
integral given for example by equation (4) in Fromang et 
al. (2002). By applying the method of steepest descent, 
this dimensionless integral can be approximated by an 
analytic asymptotic formula. However, we evaluated the 
integral numerically and checked it against the value ob- 
tained by the asymptotic approximation. 4 

We calculated the X-ray optical depth rx along the 
line of sight between the X-ray source and the point in 
question (as did Fromang et al. 2002). 5 

For a given effective X-ray ionisation rate £cff we must 
specify the fraction of H + and ions produced due to 
X-ray ionisation of H 2 . We adopted the values given by 
Cravens & Dalgarno (1978) taken from calculations of the 
efficiencies that cosmic rays with energy of 1 MeV, travers- 
ing a neutral gas of molecular hydrogen, produce these 
particular ions (see tableEJ) ■ Contributions from secondary 
ionisation were included too. Cravens & Dalgarno found 
that at proton energies of 1 MeV for example, the total 
ionisation rate (which is the sum of primary and secondary 
rates) is larger than the primary rate by a factor of 1.44. 

Since our kinetic models do not include photoreactions 
driven by incident FUV flux, we have simplified the ter- 
minology used in certain sections of this paper. Instead of 
"effective X-ray ionisation rate" we often use the phrase 
"photoionisation" . 



5. Results 

The majority of our calculations used an a-disk model 
with a = 10~ 2 and M = 10 -7 M©yr -1 . We have evolved 
the disk chemistry using the kinetic models described in 
section [31 with the primary purpose of examining the size 
of the "active zone" (the meaning of this phrase is ex- 
plained below) . One primary aim of this work is to com- 
pare and understand the differences in the size of the "ac- 
tive" zones obtained with each of these chemical networks, 
under physical conditions that are otherwise the same. 
While a number of studies of the ionisation fraction in 
protostellar disk models have appeared in the literature, 
they have each used different kinetic and/or disk models, 
so that direct comparison between them has not been pos- 
sible. 

As outlined in the introduction, the electron fraction 
is of great importance in protoplanetary disks as it deter- 
mines how well the gas is coupled to the magnetic field, 
and hence whether MHD turbulence can be sustained by 



4 Note the typo in equation (5) of Fromang et al. (2002), 
with b instead of —b. 

5 Though this particular matter is not explicitly described in 
Semenov et al. (2004), the tabulated data indicate they calcu- 
lated rx along the line of sight too. 




R [AU] 

Fig. 7. - The equilibrium electron concentration XR e [e~] 
for the transition layer separating the active from the dead 
zone vs the radial position R. A magnetic Reynolds num- 
ber of 100 is assumed. Disk regions inside R < .5 AU are 
entirely active. The disk parameters are a = 10~ 2 and 
M = 10~ 7 Mgyr" 1 . 

the MRI. The important discriminant is the magnetic 
Reynolds number Re m , which we define as 



Re m = 



(4) 



where H is the disk semi-thickness, c s is the sound speed, 
and ^ m is the magnetic diffusivity (not to be confused with 
the mean molecular weight). Numerical simulations (e.g. 
Fleming, Stone & Hawley 2000) indicate that a critical 
value of the magnetic Reynolds number, Re™ , exists such 
that non linear MHD turbulence cannot be sustained if 
Re m falls below -Re™ 1 '. The value of -Re™* depends on the 
field geometry. For a net flux vertical field -Re™' ~ 10 2 , 
and for a zero net flux field Re™ ~ 10 4 . A recent numer- 
ical study of MRI driven turbulence including the Hall 
effect indicates that -Re"/' is not particularly sensitive to 
this effect (Sano & Stone 2002). The magnetic diffusivity 
/j, m expressed in c.g.s. units is (e.g. Blaes & Balbus 1994) 



Mb 



234 



T i/2_ 



(5) 



Given the electron fraction x[e~] from the kinetic models, 
and physical conditions within the disk, we are thus able 
to calculate the spatial variation of Re m inside the disk 
models. 

In this study we take the value of -Re™ 1 * = 100. For val- 
ues of Re m less than this we assume that the disk is mag- 
netically decoupled and is "dead" with respect to sustain- 
ing MHD turbulence. For Re m > 100 we assume that the 
disk is turbulent. Following previous authors (e.g. Gammie 
1996) we use the terms "dead" and "active" zones. We 
refer to the region where Re m — 100 separating the ac- 
tive and dead zones as the "transition zone" . Values for 
the ionisation fraction x[e~] along the transition zone are 




Fig. 8. Column densities of the whole disk (solid line) and 
of the active zones (dashed and dotted lines) - referring to 
magnetic Reynolds numbers greater than 100 - for values 
x M = 0, 10" 13 ,10" 12 . The dotted line refers to x K = 
(i.e. neglecting the contribution due to the thermal ioni- 
sation of potassium). The disk parameters are a = 10 -2 
and M = 10~ 7 Moyr" 1 . 



Fig. 9. - model l/model4 - Column densities of the whole 
disk (solid line) and of the active zones (dashed and dotted 
lines) - referring to magnetic Reynolds numbers greater 
than 100 - for different values ^M(g)- While the dashed 
lines refer to a value x gT = (gas-phase ch.), the dot- 
ted lines refer to a value x gl . — 10~ 12 (gas-grain ch.). 
Especially, for a grain free environment no dead zones are 
observed above values xm > 10~ n . The disk parameters 
are a = 10 -2 and M = 10~ 7 M Q yr _1 . 



given in figure showing the values of x[e~] required to 
satisfy Re m — 100 at each radial location in the disk. 

When presenting our results we typically plot the col- 
umn density of the active zone and of the total disk mass 
as a function of cylindrical radius in the disk, rather than 
presenting contour plots of the free-electron abundance 
or magnetic Reynolds number. This is to aid comparison 
with results presented by Fromang et al. (2002) , who also 
presented their results in this way, and also because it is 
then clear what fraction of the disk mass is active at a 
particular radial location in the disk. 

We solved the kinetic equations by integrating over a 
time interval of 100,000 yrs. Hence, the ionisation fraction 
x[e~] is a function of time t, and in principle so is the loca- 
tion of the transition zone. However, for all kinetic models 
the change in the vertical location of the transition zone 
at all cylindrical radii in the computational domain was 
below the grid resolution for t > 10, 000 yrs. We find that 
after 100,000 years, the networks modell, model2, model4 
and model5 have reached chemical equilibrium. The more 
complex models, however, (model3, model6 and model7), 
do not reach strict chemical equilibrium in this time, al- 
though the position of the transition zone remains ap- 
proximately constant. A similar situation has been found, 
for example, by Semenov et al. (2004), who also consid- 
ered complex chemical models. These authors examined 
the evolution for up to 10 6 years, and found that strict 
chemical equilibrium was still not achieved in optically 
thick regions. 



5.1. Test case 

In order to test the fidelity of our disk model and X-ray 
ionisation scheme we have recomputed one of the models 
presented in Fromang et al. (2002), who used the kinetic 
model of Oppenheimer & Dalgarno (1974) to compute the 
ionisation fraction in a-disk models. Equation JQl shows 
that the equilibrium free-electron abundance Xoo[e _ ] de- 
pends on the effective X-ray ionisation rate Ceff(?", z), the 
number density of the neutral particles n, and the gas 
temperature T. We calculated the electron fraction using 
the same parameter values used by Fromang et al. (2002), 
and obtained the spatial distribution of Re m using equa- 
tions l@J and JSJ . We plot the column densities of the ac- 
tive zone obtained by our model in figure [8] which should 
be compared with figures (4) and (5) in Fromang et al. 
(2002). The agreement is excellent. 

For high temperature (and density) regions, we have 
considered the contributions of ground state ionized atoms 
under the assumption of local thermodynamic equilib- 
rium. Despite their low abundances, potassium (K) and 
sodium (Na) provide major contributions the ionisation 
fraction due to the small values of the first ionisation po- 
tentials, 5.14 and 4.34 eV, respectively. 

Applying equation (1) of Fromang et al. (2002), we 
considered the contribution of thermal ionisation of potas- 
sium to the electron abundance. Neglecting this particular 
contribution leads to smaller values in the column densi- 
ties of the active zone in the very inner regions only as 
shown in figure |H1 Including it ensures that the whole disk 
within about 0.4 AU is active. 



12 



M. Ilgner, R.P. Nelson: Ionisation fraction in disks 




Fig. 10. - model2 - Column densities of the whole disk 
(solid line) and of the active zones (dashed lines) - re- 
ferring to magnetic Reynolds numbers greater than 100 - 
for values XMg = 0, 10~ 12 , 10" 10 . The disk parameters are 
a = 10" 2 and M = 10" 7 Moyr" 1 . 

5.2. Gas-phase chemistry 
modell 

In section 15.11 the electron fraction predicted by the 
Oppenheimer & Dalgarno (1974) kinetic model was ob- 
tained using an analytic approximation, giving rise to the 
results in figure |H1 We have also solved the full kinetic 
equations of modell, and can compare the exact equi- 
librium value of x[e~] with the equilibrium value a;oo[e _ ] 
given by equation (Q. 

We find that the solution for Xoo[e _ ] obtained from 
equation Q is very accurate if the exact equilibrium val- 
ues for neutral metal atoms, x M , is used in equation ljT|). 
Comparing our results in figure El to those presented in 
figure [S] and Fromang et al. (2002), we obtain agree- 
ment when the metal abundance x M = 0. We find sig- 
nificant differences, however, when the metal abundance 
is x M — 1CP 13 or 10~ 12 . The differences can be traced to 
the assumption, used to obtain figure |H1 that the equilib- 
rium concentration of neutral metals, Xoo[M], is constant 
and equal to the total fractional abundance of metals. In 
fact the full solution to the kinetic equations predicts that 
x[M] can vary by up to two orders of magnitude, espe- 
cially within R < 2 AU, where the effects of the metals 
in increasing the vertical size of the active zone is corre- 
spondingly diminished. 

As already discussed at some length by Fromang et al. 
(2002), the addition of small quantities of metals (much 
below solar abundance) to the reaction scheme causes 
the predicted size of the active zone to increase. This is 
because the charge-transfer reaction between the metal 
and molecular ion reduces the concentration of molecular 
ions which recombine with electrons rapidly. The recom- 
bination rate between metal ions and electrons is much 
smaller, leading to a higher equilibrium electron fraction. 



Fig. 11. - model3/model7 - Column densities of the whole 
disk (solid line) and of the active zones (dashed and dotted 
lines) - referring to magnetic Reynolds numbers greater 
than 100 - for values x Mg = 0, 10~ 12 , 1CT 10 , 1CT 8 (x Fc = 
0). The dashed lines refer to a value x gr = (gas-phase 
ch.), while the dotted lines refer to x gr — 10~ 12 (gas- 
grain ch.). The disk parameters are a = 1CP 2 and M = 
10~ 7 Moyr- 1 . 



Comparing the rate coefficient for molecular ion-electron 
recombination (a) with that for metal ion-electron recom- 
bination (7) illustrates this point (see table 

Moving out from the centre of the disk, the features 
observed in the size of the active zone can be explained 
as follows. Interior to R < 0.5 AU the disk is completely 
active due to thermal ionisation of potassium. For values 
of x M = 10~ 13 and 10" 12 the region between 0.5 < R < 2 
AU maintains a deep dead zone because the metal concen- 
tration is too small to prevent Hjj" being the dominant ion, 
whose continued presence maintains a low electron frac- 
tion. Exterior to 2 AU a metal concentration significantly 
enlarges the active zone. Primarily this is because a lower 
ionisation degree is required at larger radius to obtain an 
active layer with Re m > 100, as shown by figure A 
relatively low abundance of metals is required to achieve 
this threshold at radii beyond 2 AU. Note that for metal 
concentrations x M > 10" 11 the whole of our disk model is 
active. 

mode!2 

This model takes the elements hydrogen and magnesium, 
and extracts all the relevant species and reactions from 
the UMIST database. Solving the resulting reaction set 
for mode 12 produces an active zone with much larger col- 
umn density than obtained using modell. This is illus- 
trated by comparing figures [5] and ^| It is clear that this 
effect has nothing to do with the presence of metals, as the 
active zone is large in their absence. Instead the hydrogen 
chemistry contained in the UMIST database leads to a 
dominant molecular ion H + whose recombination time is 
much longer than that (for Hj) in the Oppenheimer & 
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Fig. 12. - mode 16 - Column densities of the whole disk 
(solid line) and of the active zones (dashed and dotted 
lines) - referring to magnetic Reynolds numbers greater 
than 100 - for values x Mg = 0, 10~ 12 , lO" 10 , 10~ 8 . The 
dashed lines refer to a value x gI = (gas-phase ch.), while 
the dotted lines refer to x gI = 10 -12 (gas-grain ch.). The 
disk parameters are a — 10~ 2 and M = 10 -7 M^yr^ 1 . 



Fig. 13. - model5 - Column densities of the whole disk 
(solid line) and of the active zones (dashed and dotted 
lines) - referring to magnetic Reynolds numbers greater 
than 100 - for different values of i M ( g j . While the dashed 
lines refer to x gI = (gas-phase ch.), the dotted lines refer 
to x gI — 10~ 12 (gas-grain ch.). The disk parameters are 
a = 1Q- 2 and M = 10~ 7 Movr -1 . 



Dalgarno model. The active zone predicted by model2 is 
also much larger than predicted by a much more complex 
gas-phase chemistry (model3) so we neglect it from fur- 
ther discussion as it seems to be unreliable. 

mode!3 

This model was constructed by extracting all species and 
reactions from table IA.1I and the UMIST database con- 
taining the elements H, He, C, O, N, S, Si, Mg, Fe. The 
results obtained when solving the mode 13 reaction set are 
shown in figure 1111 Comparing this with figure El it is 
apparent that mode 13 predicts smaller active zones than 
modell, both with and without metals. The metal atoms 
included in model3 were magnesium and iron, however we 
find essentially no difference in the ionisation structure of 
the disk when iron is neglected. 

Apart from the chemiionisation reaction discussed in 
section 01 free electrons can only be created by the pho- 
toionisation of H2 and He. The inclusion of He thus in- 
creases the ionisation rate compared to that in modell. 
Removal of free electrons occurs primarily through dis- 
sociative recombination reactions with ions. Most of these 
ions are formed through charge transfer-reactions that are 
initiated by the photoionisation of hydrogen and helium. 
This leads to the formation of about 100 ionized gas-phase 
species in model3 which may recombine with electrons. 

For XMg = 0, the only ion in modell is H^. In model3 
the dominant molecular ions are HCNH + , H 4 C 2 N + , 
CNC + and NO + . Although their concentration values are 
below the value of a;[H^"] obtained using modell, they de- 
stroy e~ more efficiently than H^. 

For XMg 7^ 0, Mg + becomes by far the most dominant 
ion in the active zone as XMg is increased. However, due 



to the much larger rate coefficients for dissociative recom- 
bination reactions involving ions such as NHj, NO + and 
HCNH + , these more complex species remain the domi- 
nant species that recombine with free electrons. Hence, 
Mg + ions change rather than dominate the destruction 
rate for free electrons. This remains the case even for rel- 
atively large metal abundances up to XMg = 10~ 8 , such 
that a dead zone remains for this value. 

Motivated by the results of Sano et al. (2000) and those 
obtained using model3 which show the nitrogen bearing 
species to be important, we studied the effect of nitrogen 
on the ionisation fraction by implementing model6, which 
does not include nitrogen bearing species. We considered 
the gas-phase kinetics of model6 by setting x gI — 0. The 
molecular ions which dominate the destruction of e _ now 
differ from those in mode 13. Dissociative recombination re- 
actions with ions such as H 3 + , HCO + (for XMg = 0) and 
Mg+, H3O+, C 3 H+ and CH 3 CO+ (for x Mg ^ 0) dominate 
the destruction process of e~. However, we find that the 
location of the transition zone is only significantly modi- 
fied in the regions interior to B = 2 AU, as may be ob- 
served when comparing figures El an d El The neglect 
of nitrogen bearing species increases the size of the ac- 
tive zone, indicating that these species are important in 
determining the ionisation balance in protostellar disks. 
However, it is also clear that the remaining non-nitrogen 
bearing molecular ions are equally important, and lead to 
the formation of a smaller active zone than predicted by 
the simple Oppenheimer & Dalgarno (1974) model. Even 
a metal abundance of x Mg = 10~ 8 is insufficient to make 
the disk fully active due to the continued dominance of 
the molecular ions. 

We conclude that the additional molecular ions gener- 
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ated in more complex reaction networks generally lead to 
a smaller fractional abundance of free electrons than the 
simple Oppenheimer & Dalgarno (1974) reaction scheme. 
Hence a larger dead zone occurs within protostellar disks 
than is predicted using the simple scheme in conjunction 
with the rate coefficients listed in tabled However, in sec- 
tion|Hlwe show that by making modest changes to the rate 
coefficients in the Oppenheimer & Dalgarno model good 
agreement may be obtained with the full UMIST database 
model for certain metal abundances. 

5.3. Gas-grain chemistry 

We now present results from kinetic models that in- 
clude gas-phase reactions, "mantle chemistry" , and "grain 
chemistry" (see section 13.1.31 for a definition of these 
terms). We consider grains of a fixed size (r gr = 1CP 5 cm) 
and with the possibility of having five different charges: 
gr, gr^, and gr 2± . Our motivation for examining these 
models is to first address the question of how the vari- 
ous gas-phase models change and differ from one another 
when grains are introduced. Second, to examine the ef- 
fects of depleting the concentration of small grains on the 
chemistry as a crude means of modelling the effects of 
grain growth, and third to define under which conditions 
gas-grain chemistry allows the existence of a significant 
active zone. 

model4/model5 

The first reaction scheme including gas-grain chemistry 
used to examine the issue of the ionisation fraction in pro- 
tostellar disks was introduced by Sano et al. (2000). Their 
model is equivalent to our model5. We begin our discus- 
sion by first examining the gas-phase chemistry produced 
by model5 (i.e. with x gI = 0). The results are shown 
in figure 1131 Comparing figure 1131 with figure |H1 we find 
that model5 is in close agreement with model 1 (i.e. the 
Oppenheimer & Dalgarno model for which x[e~] is cal- 
culated from equation We note that the introduction 
of grain chemistry dramatically decreases the size of the 
active zone, and gives agreement with model4. 

We now discuss the results obtained with model4. This 
model is an extension of the Oppenheimer & Dalgarno 
(1974) reaction scheme (modell) which includes gas-grain 
interactions. The original model included a generalised 
molecule and a generalised metal. As we need to spec- 
ify the evaporation temperature from the grains for these 
components, we use the species H2 and Mg. 

The results obtained with model4, but with magne- 
sium not being included, are presented in figure 1141 (we 
discuss the results with magnesium included below). In 
the limit of x gI = this model reduces to modell. The in- 
troduction of grains with concentration x gY = 10~ 12 leads 
to a dramatic reduction in the size of the active zone, as 
the grains are very efficient at sweeping up most of the 
free electrons. This result is in basic agreement with ex- 
pectations. 

We now consider the effect of reducing the grain con- 



centration. Figure ITU shows the effect of including grains 
with concentrations x gI — 10 -14 and 10 -16 . A grain con- 
centration of a;g r = 10~ 16 is required for the gas-grain 
chemistry to generate an active zone that is the same size 
as the pure gas-phase chemistry. This result is interest- 
ing, and at first sight rather puzzling, since figure \7\ shows 
that the electron fraction required for the disk to be ac- 
tive is between x[e~] = 3 x 10~ 13 and 4 x 10~ 12 . If the 
grains simply act as a means of sweeping up the free elec- 
trons, how can a very small concentration of grains such 
= 10 14 significantly affect the position of the tran- 
sition between active and dead zones? 

In the simple model4 without magnesium, the grain 
particles are involved in the main destruction paths for 
H^~ and e~. While free electrons are mostly destroyed by 
neutral grains through the reaction 

c~ + gr — ► gr~ , 

H^~ is destroyed by recombination with singly negatively 
charged grain particles 

H+ + gr" — > H 2 + gr. 

These two processes form a loop that recycles neutral 
grain particles, ensuring that they always participate ac- 
tively in the chemical evolution. The equilibrium elec- 
tron fraction ^[e - ] is determined by a balance between 
production and destruction processes. Apparently only a 
small number of active neutral grains are required to mod- 
ify the equilibrium ionisation structure in a disk due to the 
high rate with which electrons accrete onto the grains. 

In order to reduce the destruction of e~ by neutral 
grains one may limit the neutral grain concentration x[gr] 
by: 

(i) . reducing x gr - the total grain concentration 

(ii) . limit the formation of neutral grains by breaking 

the cycle described above. 
Satisfying option (i) requires the value x gI = 10~ 16 al- 
ready discussed. By contrast, we find that the transition 
from gas-grain dominated chemistry to gas-phase chem- 
istry can be achieved for values x gY ~ 10 -14 if the cycle 
is short-circuited by setting the rate coefficient of the re- 
action H^ + gr - — ► H 2 + gr to zero (i.e by preventing 
the neutralisation of negative grains). Preventing the con- 
version of gr 2- to gr - by a similar reaction allows a grain 
concentration of x gr ~ 10~ 13 to reproduce the gas-phase 
chemistry. This is because all the neutral grains become 
and remain negatively charged, after which they no longer 
participate in the chemical evolution. 

We now discuss the results of model4 with magne- 
sium being included with concentration iMg = 10 -12 . 
These calculations are presented in figure ITol The grains 
are very efficient at sweeping up magnesium, and ther- 
mal desorption from the grain surfaces is rather inefficient. 
The result is that the active zone with grain concentration 
Xgr = 10~ 12 and magnesium concentration iM g = 10~ 12 
is almost indistinguishable from the case with iM g = 
(see figure El . 




Fig. 14. - model4 - Column densities of the whole disk 
(solid line) and of the active zones (dashed and dotted 
lines) - referring to magnetic Reynolds numbers greater 
than 100 - for different values of x gv (while XMg — is as- 
sumed). While the dashed line refers to x gI = (gas-phase 
ch.), the dotted lines refer to x gT ^ (gas-grain ch.). The 
disk parameters are a — 10 -2 and M = 10 -7 Mqjt^ 1 . 

A series of runs were performed to examine the effect of 
reducing the grain concentration x gr . In particular we are 
interested in defining the grain concentration required to 
reproduce the pure gas-phase chemistry of model 1 with 
iMg = 10~ 12 . The results of these calculations are pre- 
sented in figure Reducing x gI to 10 -16 results in an 
active zone whose size is the same as that produced by 
the pure gas-phase chemistry (modell) but with XMg = 0. 
In order to obtain an active zone that is the same size as 
that obtained with the gas-phase chemistry (modell) with 
XMg = 10~ 12 we had to reduce the grain concentration 
by more than eight orders of magnitude to x gT = 10 -21 - 
10~ 20 . The reasons for this are similar to those already 
discussed for mode 14 with iMg = 0. The equilibrium abun- 
dance of free electrons is determined by a balance between 
production and destruction processes. Production occurs 
through photoionisation of H2 only. Destruction is dom- 
inated by the reaction gr + e~ — * g r_ - Only a very 
small number of neutral grains are required to reduce the 
equilibrium electron fraction below that obtained by the 
gas-phase chemistry. When Mg is included in the model, 
the gas-phase chemistry produces a larger electron frac- 
tion than when XMg = 0. Obtaining this enhanced elec- 
tron fraction means that the required depletion of grains 
is greater when XMg 7^ 0. 

Neutral grains are recycled in this case by the reac- 
tion gr~ + Mg + — > gr + Mg as Mg + is the dom- 
inant gas-phase ion. By switching off this reaction (and 
the counterpart for gr 2- + Mg + ) in model4, and all reac- 
tions involving the adsorption of magnesium onto grains, 
we obtained an active zone whose size was the same as 
that obtained by the pure gas-phase chemistry (modell), 
but now with x gl - = 10~ 16 rather than 10~ 21 . 



Fig. 15. - model4 - Column densities of the whole disk 
(solid line) and of the active zones (dashed and dotted 
lines) - referring to magnetic Reynolds numbers greater 
than 100 - for different values x gr (ijjg = 10~ 12 is as- 
sumed). While the dashed line refers to a value x gl . = 
(gas-phase ch.), the dotted lines refer to x gl . ^ (gas- 
grain ch.). The disk parameters are a = 10~ 2 and M = 
1Q- 7 M©yr _1 . 



Simply switching off adsorption of magnesium onto 
grains, but continuing to allow Mg + to neutralise nega- 
tively charged grains gr~ resulted in a value of x gr = 10 -20 
being required to reproduced the gas-phase chemistry. 
Thus the production of a sizeable active zone in the pres- 
ence of small grains clearly requires either a highly di- 
minished population of these small grains, or a means of 
preventing negatively charged grains from being continu- 
ously neutralised by positive ions. 

model6 / model7 

We now discuss the results obtained with model6 and 
model7. We remind the reader that model6 was con- 
structed by examining the elements contained in the Sano 
et al. (2000) model (equivalent to our model5), and ex- 
tracting all species from table IA.1I All reactions of the 
UMIST database are considered that involve these species, 
mo del 7 was constructed similarly by specifying a larger set 
of elements (H, He, C, N. O, S, Si, Mg, Fe). 

As a general result we find that the location of the 
transition zone is very similar when comparing model4, 
model6 and model7 when we use a standard grain concen- 
tration of x gY — 10~ 12 . This shows that the grain chem- 
istry is dominant in determining the ionisation fraction. 
Given this fact we do not discuss model6 any further. 

We first consider model7 when magnesium is neglected 
from the model. Results for this case are shown in figure lTf)l 
for a variety of chosen grain concentrations. We will dis- 
cuss cases with xyig 7^ later on in this section. 

For grain concentration a; gr = 10 -12 the dead zone is 
significantly larger than obtained with the pure gas-phase 
equivalent mode 13. Decreasing the grain concentration to 
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Fig. 16. - model7 - Column densities of the whole disk 
(solid line) and of the active zones (dashed and dotted 
lines) - referring to magnetic Reynolds numbers greater 
than 100 - for different values of x gr (while x~f c = xyig = 
is assumed). While the dashed line refers to x gI = 
(gas-phase ch.), the dotted lines refer to x gI ^ (gas- 
grain ch.). The disk parameters are a — 10~ 2 and M = 
10~ 7 Movr" 1 . 
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Fig. 17. - model7 - Column densities of the whole disk 
(solid line) and of the active zones (dashed and dotted 
lines) - referring to magnetic Reynolds numbers greater 
than 100 - for different values of x gI (while xf c — 0, iM g = 
10 -8 is assumed). While the dashed line refers to x gr — 
(gas-phase ch.), the dotted lines refer to x gI ^ (gas- 
grain ch.). The disk parameters are a = 10~ 2 and M = 
10~ 7 Moyr- 1 . 



x gI = 10~ 16 resulted in an active zone which is the same 
size as obtained by the pure gas-phase model (model3). 
The reasons for this were discussed at length when describ- 
ing the results from model4, and we do not repeat them 
here. We have run numerical experiments in which singly 
negative-charged grains are prevented from neutralising 
via recombination reactions with positive ions. In these 
models a grain concentration of x gr ~ 10~ 13 was required 
to reproduce the size of the active zone obtained in the 
equivalent pure gas-phase model (model3). 

We now discuss results obtained using mode 17 with 
%Mg = 10 -8 . These are presented in figure IITI For a grain 
concentration of x gI — 10~ 12 we find that the active zone 
is almost identical to that obtained when i^g = 0. This 
is because the grains are very effective at sweeping up the 
metal atoms and free electrons. As with model4, a reduc- 
tion of the grain concentration to x gI = 10~ 16 leads to an 
active zone whose size is very similar to that obtained us- 
ing a pure gas-phase chemistry without metals (model3). 
A further reduction is grain concentration to x gv = 10 -20 
is required to obtain an active zone of the same size as 
obtained by pure gas-phase chemistry with magnesium 
abundance XMg = 10~ 8 . Once again the reasons are the 
same as described for model4, and we do not repeat them 
here. However, it is worth re-emphasising the fact that 
the active zone obtained by the pure gas-phase chemistry 
of model3 can be obtained for a smaller grain depletion 
factor (:Eg r = 10~ 16 ) if singly-negative charged grains are 
prevented from forming neutral grains via recombination 
reactions with positive ions. 



6. Discussion 

We now discuss issues relating to potential omissions 
from and extensions to the kinetic models presented in 
sectional Apart from the model discussed in "Fitting the 
x[e _ ] distribution" all other models we discuss in this 
section assumed XMg = 1-09 x 10 -8 and xf = 2.74 x 10~ 9 
instead of XMg — variable and xf c = 0. 

Radionuclides 

Sano et al. (2000) and Semenov et al. (2004) included 
ionisation due to the decay of radioisotopes in their 
models. We have run models with this included, us- 
ing a conservative estimate of the ionisation rate of 
Cr w 6.9 x 10~ 23 s- 1 following Sano et al. (2000). We 
find that this is not an important source of ionisation in 
our models. 

Nonthermal desorption processes 

Najita et al. (2001) examined the potential effect of 
nonthermal desorption of grain mantles by X-ray irra- 
diation. We have investigated this effect in the simplest 
possible way, by switching off all reactions that involve 
the adsorption of species onto grain surfaces. It should 
be noted, however, that we still allow electrons to accrete 
onto grain surfaces at the usual rates. We re-ran model7 
with this modification to the reaction network. The 
results are shown in figure 1181 demonstrating clearly 
that the adsorption of molecular, atomic, and ionic 
species onto grain surfaces plays essentially no role in 
determining the size of the active zone. 

Fitting the x[e~] distribution 
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Fig. 18. - model7 - Column densities of the whole disk 
(solid line) and of the active zones (dotted and dashed 
lines) - referring to magnetic Reynolds numbers greater 
than 100. While the dashed line refers to simulation which 
takes account of all types of reactions, the dotted line 
refers to a simulation where all reactions which may form a 
mantle species were neglected. The assumed metal abun- 
dances are: a; Fo = 2.74 x 10~ 9 and x Mg = 1.09 x 10~ 8 

while x„ = 10~ 12 . The disk parameters are a = 10 -2 and 
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Fig. 19. Column densities of the whole disk (solid line) 
and of the active zones (dashed and dotted lines) - re- 
ferring to magnetic Reynolds numbers greater than 100 
- for two different kinetic models of gas-phase chemistry 
with % g = 10 _12 ,10 -10 . The dashed line refers to the 
UMIST model (model3) while the dotted line refers to the 
Oppenheimer & Dalgarno model (model 1) with modified 
rate coefficients a and 7. Compared with the values in ta- 
ble^ ot and 7 are increase by 10 and 200, respectively. The 
disk parameters are a = 10~ 2 and M = 10 -7 M Q yr _1 . 



In section we computed a variety of kinetic models for 
the primary purpose of examining the size of the active 
zone. Considering first the gas-phase chemistry, we found 
significant differences in the predictions of the simplest 
model (modell) and the most complex model (model3). 
Future research in the dynamics of protoplanetary 
disks will eventually involve a direct coupling between 
magnetohydrodynamic simulations of turbulent disks, 
and chemical networks that self-consistently calculate 
the ionisation fraction. The task of coupling these two 
approaches would be greatly simplified if a chemical 
network containing few species could be developed that 
calculates the ionisation fraction accurately. For this 
reason we are motivated to examine if simple changes 
to the rate coefficients in the Oppenheimer & Dalgarno 
(1974) network (modell) can lead to agreement with 
the predicted ionisation fraction given by the complex 
UMIST-based model (model3). 

The results of our attempts to fit model3 by modifying 
modell is shown in figure 1191 For a metal abundance 
#Mg = 10~ 12 we obtain a good fit to the size of the active 
zone by increasing the rate coefficients a and 7 by factors 
of 10 and 200 compared to the values given in table ^ 
The spatial distribution of the electron fraction obtained 
using model3 is shown in figure 1201 and that obtained 
using the modified modell is shown in figure |2"T1 Although 
there are small differences, overall these distributions 
agree well. 

The agreement is not so good when we increase the 
metal abundance to XMg = 10~ 10 . We found that simple 



modifications to a and 7 were unable to produce a good 
overall fit. One is able to fit either the inner (R < 2 AU) 
or outer region (R > 2 AU) quite well, but not both at 
the same time (figure ^| shows a case where we have 
fitted the outer region). This is because the chemistry in 
modell is dominated by the metal ion when XMg > 10~ 12 , 
whereas this is not the case for model3. 

The situation with the gas-grain chemistry is much 
more simple. For a standard abundance of submicron 
sized grains, the simple model (model4) and the complex 
model (model7) are in extremely good agreement. 

Disk Model 

In addition to computing the ionisation fraction in disk 
models with a = 10~ 2 and M = 10~ 7 M©yr -1 , we also 
computed some cases for which the disk parameters are 
a = 5x 10~ 3 and M = 10~ 8 M Q yr _1 . It may be argued 
that these are closer to the canonical values for T Tauri 
stars. The effective X-ray ionisation rate for this model 
is shown in figure The lower mass and density of this 
disk model mean that the ionisation rate is now higher. 
The mass of the disk within the computational domain is 
now about 0.0049 M Q compared with 0.0087 M Q for the 
heavier disk model. 

By applying model7 we calculated the ionisation 
fraction for x gr = and x gl - = 10" 12 . The column 
densities of the active zone are shown in figure [2H We 
simply comment that a greater percentage of the disk 
is active when the surface density is reduced, but when 
grains are absent even a relatively large abundance of 
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Fig. 20. This contour plot shows the electron distribution 
at t = 10 5 yrs for our a = 1(T 2 , M = 10~ 7 Meyr" 1 
disk model by applying the UMIST model (model3) with 
iMg = 10 -12 . The contour lines refer to values x[e~] of 
1CT 14 , 10~ 13 , 10" 12 , and 10 _n . 

metals is unable to render the disk completely active. 
When grains are present, however, it remains the case 
that only about 1 % of the matter beyond R ~ 1 AU is 
magnetically active according to our adopted criterion. 



7. Summary 

We have presented calculations of the ionisation fraction 
in ct-disk models using a number of chemical reaction net- 
works that have appeared in the literature. The primary 
aims arc: to compare the predictions of these networks for 
the ionisation degree in identical disk models; to examine 
the role of gas-grain chemistry in determining the ionisa- 
tion fraction in protoplanetary disks; to examine the level 
of grain depletion required to converge toward the pure 
gas-phase chemistry. Our main conclusions are: 

(1) . Solving equation (1) to obtain the electron fraction, 
using the approximation that the neutral metal abundance 
is constant, leads to an overestimate of the size of the ac- 
tive zone when heavy metals are included in the gas. 

(2) . In agreement with previous work (Fromang et al. 
2002), we find that the simple Oppenhcimer & Dalgarno 
(1974) reaction network predicts globally active disks 
when the elemental abundance of metals exceeds 10~ n 
and the disk parameters are M = 10~ 7 M yr _1 and 
a = 10" 2 . 

(3) . A complex gas-phase chemistry drawn from the 
UMIST data base predicts larger dead zones than the 
simple Oppenhcimer & Dalgarno (1974) model. Whereas 
disks can be rendered fully active by the addition of gas- 
phase heavy metals to the simple model, the complex 
model predicts extensive dead zones even for XMg = 10~ 8 . 

(4) . A procedure for fitting the results of the complex gas- 
phase model, using a modification of the Oppenheimer & 
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Fig. 21. This contour plot shows the electron distribution 
at t = 10 5 yrs for our a = 10~ 2 , M = 10~ 7 Meyr^ 1 disk 
model by applying the Oppenhcimer & Dalgarno model 
(modell) with XMg — 10~ 12 and modified rate coefficients 
a and 7. Compared with the values in table ^ & an d 7 
are increase by 10 and 200, respectively. The contour lines 
refer to values x[e~] of 10~ 13 , 10 -12 , and 10 -11 . 



Dalgarno scheme, can be successful for specific models, 
but breaks down when the abundance of gas-phase met- 
als is increased above XMg = 10~ 12 . 

(5) . We find that the addition of small grains with con- 
centration x gr = 10~ 12 leads to a dramatic increase in the 
size of the dead zone for all gas-grain chemical models 
considered. 

(6) . In contrast to the gas-phase chemistry, there is very 
good agreement in the predicted size and structure of the 
active zone between the different gas-grain chemical net- 
works. This is because the grains play a dominant role in 
all gas-grain schemes. 

(7) . A grain depletion factor of ~ 10~ 4 is required to re- 
produce the metal-free gas phase chemistry. 

(8) . A grain depletion factor of ~ 10~ 8 is required to re- 
produce the gas phase chemistry including metals. 
Taken at face value, conclusions (7) and (8) suggest that 
efficient growth and settling of small grains will be re- 
quired if gas-phase chemical models are to be used to 
calculate the ionisation structure of protoplanetary disks. 
Furthermore, this requirement must also be met before 
a significant mass fraction of a protoplanetary disk can 
sustain MHD turbulence. In a recent study, Dullcmond 
& Dominik (2005) calculated the growth and settling of 
grains in laminar and turbulent disks, demonstrating that 
theory predicts the rapid removal of small, submicron 
sized grains. Depletion factors of < 10~ 6 were obtained 
within ~ 10 4 yr, with the consequences for observation 
and theory being: the disk is rendered optically thin to 
the UV (which has potential implications for the chem- 
istry and ionisation fraction); the 10 /im silicate feature 
is reduced or removed from the spectral energy distribu- 
tion; the mid-IR emission undergoes rapid decline. The 
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Fig. 22. The effective X-ray ionisation rate £efF per hy- 
drogen nucleus. The disk parameters are a = 5 x 10 -3 
and M — 10~ 8 M Q yr _1 . The contour lines refer to values 
of Ccff: 10" 19 and 1(T 21 s _1 . 

latter two consequences are at odds with T Tauri disk 
observations (Dullemond & Dominik 2005), which main- 
tain submicron grain populations for longer than 10 6 yr. 
A plausible scenario is that an initially quiescent disk un- 
dergoes grain growth and settling, leading to the develop- 
ment of MHD turbulence whose strength is modulated by 
the gas-phase electron fraction. Turbulent motions may 
contribute to the fragmentation of solids, as well as grain 
growth, so that a population of small grains is maintained 
in the disk. Feedback between grain size distribution and 
the strength of the turbulence may lead to a quasi-steady 
state in which the rate of growth and fragmentation of 
solids is controlled by the turbulence, whose strength is 
in turn controlled by the grain size distribution. Further 
work is clearly required to test the validity of such a sce- 
nario. 

In this paper we have concentrated on calculating quasi 
steady-state chemical and ionisation profiles in simple 
disk models, ignoring potentially important effects such as 
turbulent mixing, diffusion, and a time dependent X-ray 
flux. We have already constructed models including some 
of these effects, based on the reaction networks that we 
have described in this paper. These more complex models 
are described in Ilgner & Nelson (2005). 
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Appendix A: 

A.l. List of gaseous components and grain particles 

We considered the following elements: H, He, C, N, O, Mg, 
Si, S, Fe: 

A. 2. Kinetic model of the adsorbed components 

Collisions of gaseous components with grain particles may 
lead to adsorption of the gaseous components onto the 
surfaces of the grain particles. According to the basic as- 
sumption of kinetics, the adsorbed counterpart gX of the 
gaseous species X is considered as an individual compo- 
nent. In addition, the adsorbed species gX is characterized 
by the grain charge itself. Restricting to neutral and sim- 
ply charged grains, the adsorbed species X(gr), X(gr _ ), 
and X(gr + ) have to be considered as individual compo- 
nents so far. However, it becomes important if processes 
between ionized gaseous species, X + , and grain particles 
are described by the kinetic model. 

Processes onto grain surfaces which involve neutral and 
ionized gaseous species are listed in tabled Note, the nota- 
tion, e.g., X(gr~) refers to an adsorbed species of a neutral 
gaseous counterpart X onto a grain with simply negative 
charge. In addition be aware that not every reaction listed 
in table |3 may be applied for each gaseous species; cf. for 
example NTTJ . 

Three different processes are considered: i) adsorption 
processes - reactions 1 to 3 - , ii) desorption processes - 
reactions 6 to 8, and iii) collisional charging processes - 
reactions 4 and 5. 

Applying the mass action kinetic deterministic model 



Table A.l. List of gaseous species including free electrons 
and grain particles with different electric excess charges. 
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of the reactions listed in table I A. 31 one can derive a sin- 
gle ordinary differential equation for each component. The 
number of ordinary differential equations can be reduced - 
which becomes important for systems with a large number 
of gaseous components - if one considers the total contri- 
butions of the mantle components only. Taking 




(A.l) 
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Table A. 2. Evaporation energies of the considered com- 
ponents adsorbed onto grain surfaces. The values are tabu- 
lated in Hasegawa & Herbst (1993), if not referred to else- 
where. Values marked by * and f are listed in Sandford & 
Allamandola (1993) and Yamamoto et al. (1983), respec- 
tively. 



with n[gX] = n[X(gr~)] + n[X(gr)] + n[X(gr+)] and n gr = 
n[gr~] + n[gr] + n[gr+]. The quantity n[. . .] refers to the 
particle number density per cm 3 . 
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Table A. 3. Mantle chemistry 



1. X + gr" 

2. X + gr 

3. X + gr+ 

4. X+ + gr" 

5. X+ + gr 

6. X(gr-) 

7. X(gr) 
8- X(gr+) 
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into account (cf subsection 3.2) one gets 

-|n[X] = fc 2 n[X+]n[gr-] + fc 4 n[gX] - 
fcm[X]n gr 

^n[X+] = - k 2 n[X+]n[&-] -k 3 n[X + ]n[&] 

^■n[gX] = fcm[X]n gr + fc 3 n[X+]n[gr] - 
fc 4 n[gX] 

^n[gr"] = - I 2 n[X+]n[gr~] 

^n[gr] = k 2 n[X + ]n[gT^} - I 3 n[X + ]n[gr] 

-^n[gr+] = fc 3 n[X+]rc[gr] , 



